Tutorials for SpaCE (Spatial Cellular Estimator)
Here, we will show how to use SpaCE to deconvolve tumor spatial transcriptomics (ST) data, visualize the deconvolution results, and explore the mechanism of cell-cell interaction during cancer progression.
Load libraries
In order to run SpaCE, first you need to load SpaCE and other libraries required into your R session.
library(SpaCE)
library(Matrix)
library(jsonlite)
library(png)
library(ggplot2)
library(cowplot)
library(GSVA)
library(infercnv)
library(Polychrome)
library(ComplexHeatmap)
library(psych)
library(ppcor)
library(ggraph)
library(tidygraph)
library(parallel)
Create SpaCE object
To read their ST data into R, user can create an SpaCE object by create.SpaCE. Basically, create.SpaCE requires two parameters seqPath and imagePath, which point the standard output folders of 10x Visium data.
The seqPath folder should include
“barcodes.tsv.gz”: spot level barcodes
“features.tsv.gz”: list of genes
“matrix.mtx.gz”: (sparse) matrix of counts
The imagePath folder should include
“tissue_positions_list.csv” : barcodes and spatial information
“tissue_lowres_image.png” : hematoxylin and eosin (H&E) image
“scalefactors_json.json” : scaling factors for adjusting the coordinates
Here, we use a breast cancer ST dataset from 10x Visium website as an example.
library(SpaCE)
# load the in-house breast cancer ST data. User also can set these two pathes to their own data.
seqPath <- system.file("extdata",'seq',package = 'SpaCE')
imagePath <- system.file("extdata",'image',package = 'SpaCE')
# load ST data to create an SpaCE object.
ST <- create.SpaCE(seqPath,imagePath)
The count matrix is stored in ST@input$counts.
#show the ST counts (gene × spot)
ST@input$counts[1:8,1:3]
# 394.480194075x431.3861343 446.88118365x150.7920777 182.747522925x404.727718725
# MIR1302-2HG 0 0 0
# FAM138A 0 0 0
# OR4F5 0 0 0
# AL627309.1 0 0 0
# AL627309.3 0 0 0
# AL627309.2 0 0 0
# AL627309.4 0 0 0
# AL732372.1 0 0 0
Show key quality control (QC) metrics
Once you have created an SpaCE object, you can use ST.metrics to compute and plot both UMI and gene counts across ST spots, respectively. Of note, the metrics matrix is stored in ST@results$metrics.
# compute the QC metrics.
ST <- ST.metrics(ST)
# plot the QC metrics.
ST.metrics.plot(ST, cols = c("lightblue", "blue", "darkblue"))
Deconvolve ST data
By two steps, ST.deconvolve would deconvolve all mixtures of ST spots into malignant, immune, and stromal cells as well as an unidentifiable component. We first estimate malignant cell faction by modeling the segmental copy number variations (CNVs) from the expression profiles of ST data. Of note, this strategy avoids the workload of building malignant references case by case. Subsequently, based on an in-house hierarchical cell lineage derived from single-cell RNA-seq data sets from diverse cancer types, we use a constrained regression model to determine stromal and immune cell fraction, and also include an unidentifiable component to calibrate cellular density variations across spatial regions.
The inferred cell-type fraction and CNV value of all ST spots are stored in ST@input$Fraction and ST@input$CNV. Of note, the latter one will be used to explore the malignant cell clone structure in the following section.
# set the path to temporary file generated during deconvolution
tempFilePath <- "/Users/rub2/Downloads/"
# deconvolve ST data
ST <- deconvolve.ST(ST,tempFilePath)
# show the ST deconvolution results (cell type × spot)
ST@results$fraction[,1:3]
# 394.480194075x431.3861343 446.88118365x150.7920777 182.747522925x404.727718725
# Malignant 1.145085e-01 6.374932e-01 9.242377e-02
# B cell 1.266962e-06 2.729358e-02 2.442029e-03
# CAF 2.810727e-01 5.567629e-02 2.113566e-01
# Dendritic 1.229707e-06 2.418508e-08 3.097959e-05
# Endothelial 8.303732e-04 8.125009e-02 7.286822e-02
# Macrophage 2.749378e-01 7.318245e-02 1.130819e-01
# Mast 3.261972e-05 1.866764e-03 6.943601e-05
# Neutrophil 2.149391e-06 1.081567e-07 9.328643e-09
# NK 5.838546e-07 3.093022e-05 2.180547e-07
# T CD4 4.777462e-06 6.065937e-03 4.760630e-02
# T CD8 8.035477e-03 4.022806e-05 4.899418e-07
# Unidentifiable 3.205725e-01 1.171004e-01 4.601200e-01
# M1 1.357387e-01 4.780340e-02 6.475622e-02
# M2 1.391991e-01 2.537905e-02 2.276754e-02
# Macrophage other 5.755968e-11 2.587107e-09 2.555817e-02
# T CD4 naive 0.000000e+00 2.018646e-03 8.440469e-03
# T helper 0.000000e+00 2.018646e-03 4.323932e-08
# Treg 0.000000e+00 2.018646e-03 3.711067e-07
# T CD4 other 1.000000e-05 1.000000e-05 3.916542e-02
Notes: CAF: cancer-associated fibroblast; NK: natural killer cell; Treg: T regulatory cell. Cell type “Macrophage other” are macrophages that do not significantly express marker genes of both macrophage M1 and M2. Cell type “T CD4 other” are T CD4 cells that do not significantly express marker genes of T CD4 naive, T helper, and Treg cells.
Visualize the ST deconvolution results
We provide two functions ST.deconv.res.violin and ST.deconv.res.scatter to visualize the ST deconvolution results. Users can summarize the cell fraction of each cell type across all ST spots by ST.deconv.res.violin and hence learn the cell composition of cancer tissue for ST profiling. For example, the current breast cancer tissue was shown to mainly consist of malignant, macrophage, CAF, endothelial, and B cells. Moreover, ST.deconv.res.scatter can be used to show the spatial distribution of cell fraction for the cell type of interest by setting the parameter cellTypes.
# summarize the deconvolution results in violin plot
ST.deconv.res.violin(ST)
# visualize the spatial distribution of cell types
ST.deconv.res.scatter(ST, cellTypes=c("Malignant","B cell","CAF") )
ST.deconv.res.scatter(ST, cellTypes=c("Macrophage","M1","M2") )
Identify the clone structure of malignant cells
To assess the spatial organization of tumor clones, we focus on the ST spots with high fractions (>0.5) of malignant cells based on the deconvolution results. Users can run ST.malignant.clone.heatmap to cluster the inferred CNV values of malignant cell spots to determine the clone structure. Of note, cutoffMalignant is used to set the cutoff of high malignant cell fractions, and nClone is the clone (cluster) number of hierarchical clustering. The estimated clone identities are stored in ST@results$clone.
# cluster the cnv values to identify malignant cell clone
ST <- ST.malignant.clone.heatmap(ST, cutoffMalignant=0.5, nClone=3)
# show the clone (cluster) identity for malignant cell spots
ST@results$clone[1:16]
# 446.88118365x150.7920777 529.35643035x231.757423425 459.133658775x414.3564315 482.549500125x367.004946825
# "3" "3" "2" "2"
#135.519800625x303.440591025 306.1386108x353.81187765 482.599005075x394.034649525 488.36633175x323.0693037
# "1" "1" "2" "2"
# 476.237619x127.07920665 517.97029185x441.2871243 158.935641975x249.331680675 488.06930205x160.866335025
# "3" "2" "1" "3"
# 347.400986625x411.1881147 282.82177935x475.519797225 312.1782147x445.0495005 411.9801939x353.61385785
# "1" "1" "1" "1" Moreover, you can obtain the spatial distribution of distinct malignant cell clones by running ST.malignant.clone.scatter.
# plot the spatial distribution of malignant cell clones
ST.malignant.clone.scatter(ST)
Explore the intercellular interaction
After decomposing cell fractions from tumor ST data, SpaCE will infer cell-cell interactions (CCIs) in tumor microenvironment based on cell co-localization and ligand-receptor co-expression analyses. We first calculate the correlation between cell fractions across ST spots to evaluate co-localization of cell-type pairs. Cell distance co-localization does not directly indicate physical interactions. Furthermore, for the co-localized cell-type pairs, we test the enrichment of up-regulated ligand and receptor gene pairs as evidence of physical interactions.
First, ST.CCI.colocalization can analyze the co-localization of cell-type pairs by calculating their Spearman correlation, and the corresponding results are stored in ST@results$colocalization. Of note, partial correlations, corrected for malignant cell fractions, were computed between non-malignant cell types.
# calculate and visualize the cell-cell co-localization
ST <- ST.CCI.colocalization(ST)
Then, we further identified cell communication evidence by ligand-receptor (L-R) interactions. We computed an L-R interactions network score for each spot as the sum of expression products between L-R pairs, normalized by the average L-R network score from 1,000 random networks with the same degrees. The p-value was calculated as the fraction of random network scores that exceeded the original score.
We implemented the aforementioned process by designing ST.CCI.LRInteraction function that can calculate and visualize the L-R interaction. The corresponding results are stored in ST@results$LRInteraction.
# calculate and visualize the L-R interaction
ST <- ST.CCI.LRInteraction(ST)
# show the L-R network score and p value for the first six ST spots.
ST@results$LRInteraction[,1:6]
Analyze the cell type pair of interest
Here, we take the cell type pair of CAF and M2 as an example.
# analyze the interaction between CAF and M2 cells
ST.CCI.cell.pair(ST.deconv.res,imagePath, c("CAF","M2"))
Tumor-immune interface
Besides the interaction significance, we found an enrichment of CAF-M2 interactions close to boundaries between tumor and immune/stromal regions. The distance between CAF-M2 and the tumor-immune border was calculated by averaging the distances between each CAF-M2 interaction spot and its nearest tumor border spot. We randomly selected the same number of spots as CAF-M2 spots from the non-malignant regions and calculated their distances to the border as the null distribution. The result showed that CAF-M2 interaction spots are significantly closer to the tumor-immune boundaries.
# compute the distance of CAF-M2 to tumor border
ST.CCI.border(ST.deconv.res,imagePath,c("CAF","M2"))
